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We study the global symmetries of SU(2) gauge theory with TV flavors of staggered fermions in the presence of 
a chemical potential. We motivate the special interest of the case N = 1 (staggered) with fermions in the adjoint 
representation of the gauge group. We present results from numerical simulations with both hybrid Monte Carlo 
and the Two-Step Multi-Bosonic algorithm. 



1. INTRODUCTION 

QCD at finite baryonic density is expected to 
display a rich phase diagram JlJ. Unfortunately 
numerical simulations at non zero density are still 
prohibitive because of the sign problem (for a re- 
view of the recent progress on that matter see 
(2)). The finite density sector can be studied nu- 
merically in a class of theories which includes any 
SU(2) gauge theory. In such cases, in fact, the 
fermionic determinant is real (and positive for 
an even number of flavors). All these theories 
can also be studied by means of Chiral Pertur- 
bation Theory (%PT) even in the finite density 
region ||. Here xPT predicts the presence of 
some Goldstone modes associated with sponta- 
neous chiral symmetry breaking which are sensi- 
tive to the chemical potential, and can therefore 
be called diquark or baryonic Goldstone modes. 
These modes drive the onset transition between 
the low and the high density regime which ap- 
pears when \x is of the order of the pion mass 
TOtt ~ y/m. The same features that ensure that 
the fermionic determinant is real are also related 
to the prediction of an early onset transition. It is 
therefore interesting to note that there is a special 
case: the SU(2) gauge theory with one flavor of 
staggered fermions in the adjoint representation 
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of the gauge group (SU(2)^- lstags g) which has 
a real fermionic determinant, but where no bary- 
onic Goldstone modes are expected. It is there- 
fore an ideal model to explore the influence of the 
sign on the observables and on the pattern of sym- 
metry breaking. We will give numerical evidence 
that this model, if constrained to the sector with 
positive determinant, still belongs to the class of 
theories with baryonic Goldstone modes. Then 
we show how the inclusion of the sign changes this 
picture, suggesting a delayed onset transition. 

As a by-product we gain important experience 
with two numerical algorithms when fi ^ 0: the 
hybrid Monte Carlo (HMC) and Two-Step Multi- 
Bosonic algorithms. 

2. SYMMETRIES 

We first present the continuum and lattice for- 
mulation of Two Color QCD. The fermionic part 
of the action for Nf flavors / in the continuum is: 

S F = df'x^ (x)(d v +igA v + fj,8 v0 )j v ip f (x) 

J (1) 

The lattice action for staggered fermions is 
Sf =J2x P (x)(D xy [U,n}+m5 xy )x P (y) 

(2) 
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where the index p runs over N flavors of staggered 
quark, and D is given by 



D *y = ^ 1 ^—( u A x ) s x, y -i>-Ul(y)6 Xt y +c 



+^(e"[/o(,)^^ 6 - e -"[/ t (y)^, 9+ o) 



(3) 



X, X are single spin component Grassmann ob- 
jects, and the phases rj^x) are defined to be 
(— l^cH t-^c-i. For brevity, we have here used 
U to denote both the fundamental (2x2 complex 
unitary) and adjoint (3x3 real orthogonal) link 
fields. 

2.1. Sign of the determinant 

As is well known, Two Color QCD has a real 
fermionic determinant also at non zero density. In 
fact the SU(2) group only has real or pseudoreal 
representations. For all such representations it is 
always possible to find a unitary operator T such 
that KT commutes with the fermionic matrix M, 
where K is the complex conjugate operator. If 
such an operator exists one can show d|| that all 
the eigenvalues of M appear in complex conjugate 
pairs and thus the determinant of M must be 
real. Although for the theories we are considering 
here det M is real, it need not be positive definite. 
Indeed, it turns out that some of these theories 
have a sign problem. It can be shown that a 
sufficient condition for a positive definite action is 
that there is one choice of T such that (KT) 2 = 
-1. 

For the continuum action, given by ([!]), one 
finds 

fundamental: T = C^®t 2 , (KT) 2 = 1 ; 

adjoint: T = C75 <X> I, (KT) 2 = -1 . u\ 

Here, C is the charge conjugation matrix. For 
staggered lattice fermions M is given by (H||), 
and 

fundamental: T = r 2 , (KT) 2 = -1 ; 
adjoint: T = I, (KT) 2 = 1 . 

We thus have a proof that the functional inte- 
gral measure is positive definite for continuum ad- 
joint quarks and fundamental staggered quarks. 



There is no such proof for continuum fundamen- 
tal quarks and staggered adjoi nt q uarks, and as 
we shall demonstrate in section 4.4, there are in- 



deed isolated real eigenvalues and hence a sign 
problem for the adjoint staggered model at large 
chemical potential fi. 

2.2. Symmetry breaking pattern 

In the chiral limit, the action has a \J(N)l <g) 
U( N) a symmetry, which for staggered fermions 
is manifest as independent V(N) symmetries for 
the even and odd sites. At \l = this enlarges to 
a XJ(2N) symmetry. This can be seen most easily 
by introducing new fields, 



Xe = (Xe, -XeT 2 ) X a = 

for fundamental quarks, and 

X e = (Xe,X t /) *o=fer 

\Xo 



Xo 
-T2X1 1 



(G) 



(7) 



for adjoint quarks. The action can then be writ- 
ten as 



S =l E Vu(x)(L + (x)-L_(x)) 



(8) 



where 



L+=X e (x) 



p^.o U u (x)X (x+u) 



L =.Y,(»(' e ^ Sv0 )ut(x-v)X o (x~v) 



In the continuum, the equivalent fields are 



fundamental: 



adjoint: 



tpL 

which gives the lagrangian 
£ = i1' t (7„(A,- ( uB i ,)4< 
where 



B„ 



BSvo 



B 



1 

-1 



(9) 
(10) 

(11) 
(12) 

(13) 
(14) 
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The chiral condensate can be written in terms of 
the new fields, 



XX : 



x e 



±1 



2 X * 



±1 



T x 

2 X ° 



(15) 



where the + sign is for fundamental fermions and 
the — sign for adjoint, while T is the unitary 
operator defined in section 2.1. A nonzero chi- 



ral condensate thereby breaks down the \J(2N) 
symmetry to 0(2N) for fundamental fermions 
and Sp(2iV) for adjoint fermions, giving rise to 
N{2N +1) and N{2N - 1) Goldstone modes re- 
spectively. Of these, there will be N 2 mesonic 
states, while the remaining N(N ± 1) will be di- 
quarks. From this we see that in the case of N = 1 
adjoint fermions, and only in this case, are there 
no diquark Goldstone modes. 

For m 0, all states remain degenerate, gain- 
ing masses m n ot ^pm in accordance with stan- 
dard PCAC arguments. As the chemical poten- 
tial fj, increases from zero, a ground state contain- 
ing baryonic matter is promoted, signalled by a 
non-zero value for the baryon number density 



1 



+ e-^ t (.T-6) X (x-6)]). 



(16) 



At zero temperature n thus serves as an order 
parameter for an onset phase transition occur- 
ing at some fi separating the vacuum from a 
state containing matter. A naive energetic ar- 
gument would suggest that the onset transition 
should occur for a value of \x equal to the mass 
per baryon charge of the lightest particle carry- 
ing non-zero baryon number. For the models dis- 
cussed in the previous paragraph in which some 
of the Goldstone modes are diquark states, those 
states will be the lightest baryons in the spec- 
trum. This means that for most variants of Two 
Color QCD we expect fi a ~ m n /2, in contrast to 
the much larger value toat/3 expected in physical 
QCD. The exception is SU(2) 



jV = lstas 
Adj 



3. ALGORITHMS 

We have studied Two Color lattice QCD with 
N = 1 adjoint flavors of staggered fermions, using 



two different simulation algorithms: the hybrid 
Monte Carlo (HMC) algorithm §], and a Two- 
Step Multi-Bosonic (TSMB) algorithm §. In |] 
we described in detail how both algorithms are 
defined for the model above. Here we only re- 
mind that, in order to optimise the performances, 
in both cases we need to tune the simulation pa- 
rameters. In HMC we tuned the mean trajec- 
tory length and the size of the discretised time 
step. In TSMB we have a larger set of parame- 
ters: the interval [e, A] over which the polynomial 
approximation is performed, the degrees of the 
four polynomials (rij, i = 1, . . . , 4), and the num- 
ber of Metropolis, heatbath and overrelaxation 
iterations for the gauge and the boson fields re- 
spectively. 

The first observation is that as the ratio \ij 'y/rn 
increases some eigenvalues of the fermionic ma- 
trix approch zero and any simulation becomes 
hard in that region. By increasing \ij ^pm further, 
some eigenvalues get a negative real part. It turns 
out that the model allows configurations with a 
negative determinant (when an odd number of 
real negative eigenvalues appear). At this point 
TSMB and HMC display different behaviours. 
TSMB simulations can easily sample configura- 
tions with both signs of detM, for any combina- 
tion of the simulation parameters that we used. 
HMC simulations (at least for the combinations 
of trajectory lengths and time steps that we used, 
which are dictated by efficiency considerations) 
were never able to change the sign of det M. 

The reason for this difference must reside in 
the updating of the gauge fields, since det M is 
only defined in term of these. In the case of 
HMC the accept/reject step is performed only 
after updating the whole configuration and the 
exact action is used. In the case of TSMB the 
gauge field updating is performed by a Metropo- 
lis algorithm and the action is given by a poly- 
nomial approximation realised by means of aux- 
iliary boson fields. These conditions apparently 
allow TSMB to change the sign of det M more 
easily than HMC. The exactness of the TSMB al- 
gorithm is guaranteed by a final reweighting step, 
which is necessary since no polynomial approxi- 
mation is sufficient when dealing with arbitrarily 
small eigenvalues. 
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4. RESULTS 

4.1. Autocorrelation studies 

In order to establish the errors of our numeri- 
cal studies we analysed the autocorrelation times 
of both algorithms and expressed them in the 
common unit of measure of matrix multiplications 
(appropriately corrected with a factor that takes 
into account that TSMB spends more time in 
other kinds of operation). The number of ma- 
trix multiplications per sweep is a function of the 
parameters of the simulation S. It turns out, in 
general, that the gluonic observables have a much 
longer autocorrelation time than the fermionic 
ones, sometimes by even two orders of magni- 
tude. In the following we will consider only the 
former, in order to be conservative. At zero den- 
sity we could determine the autocorrelation time 
with sufficient precision for both algorithms. The 
results are shown in table [l]. The runs in the table 
are sufficiently long for the measurement of the 
integrated autocorrelations. The run in the first 
line has relatively more gauge update sweeps com- 
pared to the boson field updates: it has N M = 12 
Metropolis gauge sweeps per Nh = 2 heatbath 
and No = 1 overrelaxation bosonic sweeps. This 
is obviously not advantageous for the autocorre- 
lation. In the other two runs Nm — 4 which 
is substantially better. The difference between 
the second and third lines is in the number of 
auxiliary boson fields (equal to the rank of the 
first polynomial). The table shows that a low 
acceptance about 30% is somewhat better than 
the higher one near 70%. Further optimisation 
of the choice of TSMB parameters in this point 
is still possible but Nf^ q is already smaller than 
the corresponding number in HMC. 

Our longest TSMB run at (i = 0.4 is not long 
enough for an accurate determination of the in- 
tegrated plaquette autocorrelation. After 80000 
sweeps (on a 4 3 x 8 lattice) the obtained autocor- 
relation estimate is rf^ q ~ 2400 {N^ q ~ 6-10 6 ). 
The obtained result is good enough for an or- 
der of magnitude estimate of rf^ q but the real 
value may be somewhat larger. Also in the case 
of HMC, the longest run of 18000 trajectories was 
not sufficient to determine the autocorrelation. 

Since the point at fi = 0.4 appears so difficult, 



we started an analysis of autocorrelation at point 
/i = 0.36. Our longest run (40000 sweeps) is still 
too short to provide an accurate determination. 
However we can at least estimate the order of 
magnitude of T^ q ~ 500 (N?£ q ~ 8 • 10 5 ). 



4.2. Reweighting 

The final precision in the TSMB algorithm is 
achieved by reweighting the gauge configurations 
at the evaluation of expectation values. In a 
model with detM possibly negative, as considered 
here, the sign of the determinant is also taken 
into account in the reweighting. The choice of 
the order of the second polynomial has an im- 
portant effect on the reweighting. For sufficiently 
large 712 the two-step approximation can be so 
good that the effect of reweighting is negligible 
compared to the statistical errors. In fact, this 
ideal situation can be achieved in the low den- 
sity phase of our model where negative determi- 
nants practically never occur. In the high den- 
sity phase, however, negative determinants play 
an important role and it is more advantageous to 
have a non-negligible reweighting. (This is also 
required by the very small eigenvalues which can 
be treated exactly in the reweighting step.) Of 
course, a very low value of 712 is not optimal ei- 
ther because then the typical reweighting factors 
become too small and the effective statistics is 
substantially reduced. The effect of choosing the 
order of the second polynomial on the distribu- 
tion of reweighting factors is shown by fig. |l} 



4.3. Physics results from HMC 

We used three distinct quark masses on the 4 3 x 
8 lattice at (3 = 2.0, and explored values of \x up 
to and including 0.8 for m = 0.1, ix = 0.7 for 
in = 0.05, and fj, = 0.5 for m = 0.01. 

Here we measured the chiral condensate (xx) 
the baryon number density n (fl6|), the plaque- 
tte and the pion mass m^. If we define the 
rescaled variables x = 2fij 'm^o, V = (xx)/(xx)o, 
and n = m w o^/8m(xx)o (where the subscript 
denotes values at zero chemical potential), then 
xPT predicts (in the limit of small m and fj,) that 
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Table 1 

Integrated autocorrelation of the plaquette rf^ q for TSMB and HMC at {3 = 2.0, m = 0.1, /.i = 0.0 on 
a 4 3 x 8 lattice. For TSMB we use polynomial orders: n>i, ri2 — 90, 713 = 120; interval of approximation: 
[e = 0.005, A = 10.0]. For both algorithms R aC c is the acceptance rate of the correction step, and N sweep 
is the length of runs in sweeps or trajectories. N mu i t is the number of matrix multiplications per sweep. 
Nf^ q is the autocorrelation in number of matrix multiplications. 
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Figure 1. Examples of reweighting factors for 
different choices of the second polynomial, at 
(i = 2.0, to = 0.1, fi = 0.36 with m = 64. 



all data should fall on the lines 



x<l 
x>l 



1(1" 



x<l 

X>1 



(17) 

In Figs I and | we show the rescaled variables y 
and h respectively as functions of x. The data 
collapse very nicely onto a universal curve, cor- 
responding quite closely to the prediction (p7|). 
The systematic departures from the theoretical 



o m =0.10 
□ m=0.05 
Om=0.01 
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Figure 2. Chiral condensate vs. chemical poten- 
tial using the rescaled variables of eq. ([l7]) . 



curves for x < 2, downwards for the condensate 
data and upwards for the baryon density, may 
well be explicable by higher order corrections in 
%PT. Our more recent results from the regime 
x > 2, however, suggest a dispersion in the data 
from different m, and hence a possible breakdown 
of %PT at higher densities. This may be due to 
new thresholds as particles other than Goldstones 
are induced into the ground state, or even a fur- 
ther phase transition ||. 

Whilst the approximate quantitative agree- 
ment between our results and the theoretical pre- 
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Figure 3. Baryon density vs. chemical potential 
using the rescaled variables of eq. (p7|). 
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Figure 4. vs. fj,, for the three different quark 
masses. Also shown is the line = 2/x. 



dictions of || is gratifying, it also contradicts 
the symmetry-based arguments of section || that 
there are no baryonic Goldstones for N = 1 stag- 
gered flavor, and no gauge- invariant local diquark 
condensate. We believe that this is because the 
HMC simulations fail to take account of the de- 
terminant sign (or indeed even to change it) ie. 
that simulations with functional weight |detM| 
yield broadly similar results to those with weight 
det 2 M. The premature onset at \x Q = m w /2 is 
therefore a direct manifestation of the sign and/or 
ergodicity problems. 

In Fig. |^ we show results for , obtained us- 
ing a standard cosh fit to the meson correlator 
over all 8 timeslices. The fits are quite stable 
in the low density phase, but the correlators be- 
come very noisy once n > 0, resulting in a much 
reduced precision. It is significant, however, that 
our results in the dense phase are at least consis- 
tent with the xPT prediction = 2/j 

Finally we turn to the effect of the chemical 
potential on the gauge fields. Since this can only 
be communicated via fermion loops, any effect 
we see can be ascribed with certainty to dynam- 
ical fermions. Gluonic observables, however, are 
also much more prone to auto-correlations as de- 
scribed in section 4.1, particularly as the quark 



mass is reduced. Systematic changes with [i are 
therefore quite difficult to observe. In this initial 
HMC study we have only measured the average 
plaqucttc; the results are shown in Fig. |5|. The 
data for m = 0.1 show the plaquette remaining 
roughly constant for fj, < fi , before beginning 
to decrease. The points for m — 0.05, 0.01 have 
been omitted for clarity, but reveal a similar pic- 
ture. We interpret it as follows: for temperature 
T = 0, all values of /z < /i Q are physically equiv- 
alent corresponding to the same physical state, 
namely the vacuum. We only expect an effect 
on gluonic observables in the presence of matter, 
ie. for fi > n . To the extent that the results 
are constant for fi < [i we can be confident that 
our simulation has an effective T ~ 0. The de- 
crease in the plaquette for fi > /i may be due to 
the decrease in the number of virtual quark-anti- 
quark pairs which may form due to the Exclusion 
Principle — an effect known as Pauli blocking. 
This results in a decrease of screening via vacuum 
polarisation, and hence an effective renormalisa- 
tion of the gauge coupling (3 and consequent de- 
crease of the plaquette. In the large-^i limit the 
lattice should become saturated with one quark 
of each color per site, and the plaquette assume 
its quenched value 0]. 
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Figure 5. Average plaquette vs. n for to = 0.1. 



4.4. Physics results from TSMB 

Each TSMB simulation is characterised by a 
vector rii specifying the polynomial orders at each 
stage, as described in For each configuration 
generated, a reweighting factor r and the sign 
of det M must be determined. On the relatively 
small lattices considered here, we have been able 
to compute det M directly using standard numer- 
ical methods. The expectation value of an observ- 
able O is then determined by the ratio 



(O) 



(O x r x sign) 
(r x sign) 



(18) 



Here we present results from runs performed on 
a 4 3 x 8 lattice with (3 = 2.0, to = 0.1 and three 
values of u. The parameters are given in table [| 
Note that the polynomial orders required increase 
with ll. The point at \i = was chosen to enable 
the TSMB algorithm to be tested against HMC, 
since both should yield identical results. The 
fi =/= points were chosen so as to have one value 
just past the HMC onset transition, where the 
edge of the eigenvalue distribution just overlaps 
the line ReA = and a small percentage of nega- 
tive determinant configurations are expected, so 
that hopefully the sign problem is not too se- 
vere, and one value fairly deep in the high density 
phase. 

Our results for the standard observables, to- 



Table 2 

Parameters for the TSMB simulations. iV c f g is an 
estimate of the number of independent configura- 
tions, using the estimates of the autocorrelation 



time from section Hi. 1 
fj, 



(ni , n 2 ) iV c fg (r x sign) 



0.0 
0.36 

0.4 



(16,100) 
(64,500) 
(64,700) 
(100,1000) 



90 
120 
160 
128 



0.313(20) 
0.417(21) 
0.030(59) 



gether with the corresponding HMC results, are 
summarised in Table |. For TSMB at fx ^ 
we also include observables determined separately 
in each sign sector, defined by (0)± = (O x 
r )±/( r )±- We note that the agreement between 
HMC and TSMB at fi = 0.0 is good, as it should 
be. The results at ll = 0.4 have too large errors 
and we cannot derive any conclusion. The re- 
sults at fx = 0.36 show an acceptable agreement 
between HMC and the positive sector of TSMB. 
However, the results for the negative determinant 
sector are significantly different. (A similar trend 
may be seen at fi = 0.4, but here the difference is 
still within 2a, and thus not statistically signifi- 
cant.) The effect of the sign on the total average 
- especially for the baryon density — is quite 
interesting. Even if the errors are still large, we 
can see that n is already definitely non zero in the 
positive sector. The inclusion of the negative sec- 
tor has the effect of bringing back the average n 
to zero — although it should be pointed out that 
the TSMB average is also still only about 2a away 
from the HMC result. This mechanism is non- 
trivial since at that point there are very few con- 
figurations with negative determinant. For the 
chiral condensate, the effect of the sign is quali- 
tatively similar in that the inclusion of the sign 
counteracts the suppression observed in HMC — 
but we are unable at this stage to draw any quan- 
titative conclusions. 

These observations suggest that at this value of 
u the system is still in the low density phase, and 
hence a Q true > Mo hmc ■ This is in accord with 
our symmetry-based arguments that for N = 1 
flavors of adjoint staggered fermion there are no 
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Table 3 



A comparison of results between TSMB and HMC. At /j — 0.0 reweighting was not necessary. Due to 
the long autocorrelation times, the errors in the HMC results, especially for the plaquette, are probably 
underestimated. 
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baryonic Goldstones and hence no early onset. 
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5. CONCLUSIONS 

We presented our progress in the numerical 
study of SU(2)^= lstags . We have seen that when 
the model is constrained to the sector with pos- 
itive determinant it displays a good agreement 
with the expectations of xPT. The departure 
from xPT at large density is now more clear and 
deserves further investigation. We also confirmed, 
with narrower errorbars, the presence of a delayed 
onset when the full model is taken into account. 
We learnt how to perform a numerical analysis of 
this difficult model through the use of a TSMB al- 
gorithm which is well suited to sample correctly 
the sign. In future we shall also study the di- 
quark condensates suggested in || in order to un- 
derstand fully the effect of the sign in the region 
of the parameters where this is possible. 
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